Impact of Perturbations on Watersheds 



O 

(N 

a 



E. Fehr/'B D. Kadau,^ J. S. Andrade Jr./"'^ and H. J. Herrmann^^^ 

^IfB, ETH Zurich, 8093 Zurich, Switzerland 
Departamento de Fisica, Universidade Federal do Ceard, 60451-970 Fortaleza, Ceard, Brazil 

We find that watersheds in real and artificial landscapes can be strongly affected by small, local 
perturbations like landslides or tectonic motions. We observe power-law scaling behavior for both the 
distribution of areas enclosed by the original and the displaced watershed as well as the probability 
density to induce, after perturbation, a change at a given distance. Scaling exponents for real and 
artificial landscapes are determined, where in the latter case the exponents depend linearly on the 
Hurst exponent of the applied fractional Brownian noise. The obtained power-laws are shown to 
be independent on the strength of perturbation. Theoretical arguments relate our scaling laws for 
uncorrelated landscapes to properties of invasion percolation. 
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Watersheds are the lines separating adjacent drainage 
basins (catchments) and play, hence, a fundamental role 
in water management ^ , landslides 0, 0] and flood pre- 
vention 0, 13] ■ Since ancient times watersheds have been 
used to delimit boundaries and have already become is- 
sues in disputes between countries Moreover, sim- 
ilar problems also appear in other areas such as Image 
Processing and Medicine @, which shows the generality 
and importance to fully understand the subtle dynamical 
properties of watersheds. But how sensitive are water- 
sheds to slight localized modifications of the landscapes? 
Can these perturbations produce large, non-local changes 
in the watershed? Geographers and geomorphologists 
have studied the evolution of watersheds in time and 
found it to be driven by local events called stream cap- 
ture. These events can affect the biogeography Q, and 
occur due to erosion, natural damming, tectonic motion 
as well as volcanic activity [8-ll|. Recently, the associ- 
ated relevant mechanisms were investigated numerically 
and in small scale experiments [l2| . Finally, the problem 
studied here is also of interest to image processing, in 
order to circumvent segmentation failure 

In this Letter we investigate the effects of topological 
modifications like landslides or tectonic motion on the 
watershed. In fact, we show that the same type of topo- 
logical perturbation can indeed trigger non-local effects of 
any length scale, i.e., following power-laws distributions. 
For illustration, as shown in Fig. [1] we obtain after a lo- 
cal height change of less than 2 m at a location (cross) 
close to the Kashabowie Provincial Park, some kilometers 
North of the US-Canadian border, a substantial displace- 
ment in the watershed (blue), which encloses together 
with the original watershed (red) an area A ^ 3730 km^. 
Here a model is developed to provide a qualitative and 
quantitative description of this phenomenon. 

In our simulations, we use real and artificial landscapes 
in the form of Digital Elevation Maps (DEM), consisting 
of discretized elevation fields. Here we call sites to the 
discretization units, defined as the square areas with a 
size given by the DEM resolution. The watershed is the 




FIG. 1. (color online) Example of the watershed between US 
and Canada, close to the big lakes (red). Also shown is the 
resulting change in the watershed (blue) due to a perturbation 
of 2 m at a spot (cross) close to the border, near Thunder 
Bay. The watershed displacement encloses an area of about 
3730 km^. The dot marks the new outlet of the area after 
perturbation. The inset shows the same area on a larger scale. 

line dividing the entire landscape into two parts. Each 
part drains, according to the steepest descent along the 
coordinate directions, to either one of a chosen pair of op- 
posite boundaries (east-west or north-south) of the DEM. 
For the determination of this line we use an iterative ap- 
plication of an invasion percolation procedure (IP) [14| . 
For a given landscape, as shown in Fig. [TJ we initially 
determine its watershed (red line). Then a local event is 
induced by changing the height hk ^ hk + A at a. single 
site k (cross in Fig.[T]) of the DEM, where A is the pertur- 
bation strength. Since we are interested on the non-local 
features of the watershed response to local perturbations, 
we only perturb sites that are not on the watershed. This 
implies A > to induce changes in the watershed (blue 
line). The displacement of the watershed is quantified 
by measuring the area A between the original and the 
perturbed watershed. By definition, the water can only 
escape from the displacement area through a single site. 
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which we call outlet. The old outlet (o-outlet) before 
perturbation always coincides with the perturbed site k. 
After perturbation, k becomes part of the new watershed 
and the water escapes through a new outlet (n-outlet), 
which is located at the original watershed. After measur- 
ing the area A and the distance R between the old and 
new outlets, we proceed by restoring the original land- 
scape, i.e., the height at k is reset to its initial value. This 
procedure is repeated for every site k of the landscape, 
except those located at the original watershed. Initially, 
we fix A equal to = \hmax — hmin\, i-e., the height dif- 
ference between the lowest /imm and highest height h^ax 
of the landscape, which corresponds to a perturbation of 
infinite strength for the landscape under investigation. 
With this choice all possible changes are obtained within 
the DEM. In all definitions hereafter, we consider only 
those perturbations leading to a displacement of the wa- 
tershed. In what follows, we study the distribution P{A) 
of the areas A, the probability distribution P{R) of the 
Euclidean distance R between the two outlets, and the 
dependence between A and R. For this, we define the 
average area (A) and the distribution P{A\R) of areas A 
associated with an outlet distance R. 
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FIG. 2. (color online) The distribution P{A) is shown for 
various regions: Rocky Mountains, Andes and Appalachian 
(unshifted); Brazil and Europe (shifted by a factor of 10^ for 
better visibility); US-CAN, Kongo and Germany (shifted by 
10*). All data sets have a resolution of 540 m. The solid line 
shows the fit to the Andes data with a power-law of exponent 
— 1.65 ± 0.15. The inset shows (A) as function of R for the 
Rocky Mountains at resolutions of 270 m (squares), 540 m 
(circles) and 1350 m (triangles). The solid line has slope 2. 

First, we study several natural landscapes, from moun- 
tainous (e.g. Rocky Mountains) to rather flat landscapes 
(e.g. US-CAN, Kongo and GermanyV The DEM data 
was obtained from the SRTM-project [l5| , where for each 
set we used a size of 2700 km x 2700 km (except 1080 
kmxlOSO km for Germany), and a resolution of 540 m, 
defining the size of a site. Hence, the physical size of 
the 8 data sets are large enough, so that finite size ef- 
fects emerging from the DEM boundaries could not be 



detected. As shown in Fig. [2l we find the distribu- 
tion of areas to follow a power-law, P{A) ~ A~^ , with 
/3 = 1.65 ±0.15 for all landscapes. The probability distri- 
bution P{R) of outlet distances R also obeys a power-law, 
m?) ~ R-P, with p = 3.1 ± 0.3 « 2/3 (see Fig. SMI in 
[la|). and displays an upper cutoff in the range 50 km 
< R < 500 km for the studied landscapes. This cutoff 
is independent on the resolution and could be due to a 
length scale arising from tectonics. The value of p im- 
plies (A) ~ R^ , which agrees well with our data (inset of 
Fig. [2]). The distribution for a given distance R scales as 
P(A|i?) - A"" with a = 2.3±0.2 (see Fig. SMI in 




FIG. 3. (color online) Data collapse of the distribution P{R) 
for uncorrelated landscapes {H = — 1) of three different sys- 
tem sizes L = 129, 257, 513 (triangles, circles and squares, 
respectively). The line represents a power-law fit to the 
data for the largest landscape (squares) revealing an expo- 
nent p = 2.21 ± 0.01. The inset shows the distribution P{A) 
of the areas for the same system sizes. The line represents a 
power-law fit to the data with exponent /? — 1.16 ± 0.03. 

In order to understand these power-laws and the de- 
pendence between A and R, we study artificial land- 
scapes, where the local heights are generated using frac- 
tional Brownian motion (fBm) on a square lattice [it} . 
This model incorporates spatial long-range correlations 
to the system that are controlled by the Hurst exponent, 
H. We first consider the case of uncorrelated landscapes, 
which is a special case of the fBm model with to H = —1. 
In Fig. [21 we present the results obtained for several sys- 
tem sizes, using the same procedure as for the natural 
landscapes. The probability density P{R) again follows 
a power-law P{R) oc R"'', without upper cutoff as in real 
landscapes. We estimate p = 2.21 ± 0.01 using the scal- 
ing P{R) — LPf[RL], where L is the linear dimension of 
the landscape. For the distribution P{A) = L^f[AL'^], 
we obtain an excellent data collapse for /3 = 1.16 ± 0.03 
(see the inset of Fig. [3|). In the case of the distribu- 
tion P{A\R) at a fixed outlet distance, we again find a 
power-law P(A|i?) ~ A'" (see Fig. SM2 in [l^). Finite 
size scaling analysis yields an exponent a = 2.23 ± 0.03 
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independent on the value of R. Assuming that R de- 
scribes the extension of A in every direction, the relation 
p = a is reasonable. This is even well supported by 
the similarity of the obtained exponents. The area A 
was rescaled by L^, indicating that the areas are com- 
pact. Considering finite-size scaling, our data is consis- 
tent with (A) ~ i?2 (see Fig. SM2 in [13]). Furthermore, 
the compactness of the areas is supported by the mea- 
sured value /3 = 1.16±0.03, which agrees well with the 
relation /3 = p/2 « 1.11. 

In the following we show that we can match the ex- 
ponents quantitatively by tuning the Hurst exponent H 
to introduce spatial long-range correlations, as present in 
real geological systems. The exponents a, l3 and p were 
calculated for several values of H (see Fig. \^ . As shown 
in Fig. 01 we observe that both f3 and p increase with 
H. Furthermore, the relationship /? = p/2 is maintained, 
since the areas remain compact in the entire range of H 
values. Around H = —0.5, a starts to deviate from p 
and for if > we observe a to decrease. Previously, 
we had assumed R to reflect the extension of the area, 
i.e., the outlets to reside typically on opposite sides of 
the area. To check whether this is still valid, we mea- 
sured the angle 9 between the lines connecting the center 
of mass of the area with the two outlets. We observe 
the average angle to decrease as a function of H (see 
Fig. |4]). This implies that, on average, the two outlets 
approach each other with increasing H (see also the in- 
sets of Fig. m, so that R is no longer representative of 
the area extension. Finally we find good quantitative 
agreement with the exponents obtained from the natural 
landscapes, which are known to have a Hurst exponent 
inside the range 0.3 < H < 0.5 (see Ref. [l^ and refer- 
ences therein). Hence, except for the upper cutoff in R, 
our model provides an excellent quantitative description 
of the effects observed on natural landscapes. 

Next we analyze quantitatively the impact of the per- 
turbation strength A on the watershed. In Fig. [S] the 
number of perturbed sites that change the watershed 
is shown for uncorrelated, for artificial correlated (with 
H = 0.3) and for natural landscapes (Andes). In all 
three cases, N is found to increase linearly with the ap- 
plied perturbation strength, A. This indicates that 
changes on the watershed can be observed even for in- 
finitesimally small perturbations. Additionally, in both 
cases where correlations are present, N is observed to 
reach a plateau. As already stated, when A is equal to 
h^, this corresponds to the largest relevant perturbation, 
so that N{hj^) = Nmax ^ indicates that many per- 
turbations never change the watershed at all. It is clear 
that A > \hj — hi\ is needed, where hi and hj are the 
heights of the outlets of the area. Therefore, if the dis- 
tribution Po{h) of outlet heights is known, one obtains, 

7V(A) = J^^ £^ p,{h) poih - A') dh dA' . (1) 




-1 -0.5 0.5 



H 

FIG. 4. (color online) The exponents a (squares), (circles) 
and p (triangles) are shown for several values of the Hurst 
exponent H. Each point results from a similar study as done 
for the uncorrelated landscapes. The exponents for the nat- 
ural landscapes (open symbols), all corresponding to Hurst 
exponent values in the range 0.3 < H < 0.5, are consistent 
with our model. The average angle 6 (in radians) between the 
outlets from the center of mass is shown too (crosses). The 
insets depict schematic shapes of the areas and positions of 
the two outlets for small (left) and large (right) values of H. 




FIG. 5. (color online) Dependence of the number of perturbed 
sites A'' that promote changes on the watershed on the pertur- 
bation strength A applied for uncorrelated (squares), Andes 
(triangles) and fBm landscape with H = 0.3 (circles). The 
solid line corresponds to the analytic relation obtained from 
Eq. (m for uncorrelated landscapes. The inset shows the av- 
erage area (A) as a function of the distance R between the 
outlets for A/h^ = 1, 0.016 and 0.00025 (pluses, crosses, and 
squares, respectively), and L = 513. 

For landscapes with uniformly distributed heights, we 
find Poih) to be still a uniform distribution. Then we ob- 
tain from Eq. ©, A^(A) = ih^A-Ay2)2N^acc/iL''hl), 
which is in excellent agreement with our data (see Fig.[S]), 
where an approximately linear behavior can be observed 
for A < hw The observed power-laws are maintained 
for all values of A, as can be clearly seen for (A) in the 
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inset of Fig. [5l We conclude that infinitesimally small 
perturbations have qualitatively the same effect on the 
watershed as any larger perturbation strength A. 

In the case of uncorrelated landscapes (i/ = — 1), for a 
given area A, the corresponding invasion percolation clus- 
ter is obtained by starting the penetration process from 
one outlet to another, always growing along the steep- 
est descent. The area A can therefore be understood as 
the envelop of this IP-cluster. From percolation theory, 
the fractal dimension of the IP cluster is d/ = 91/48 in 
two dimensions jl^], which implies that (A) oc M^/"*^', 
where M is the number of sites (mass) of the cluster. 
This result is consistent with our simulations. The size- 
distribution P{M\R) of IP-clusters between two sites at 
a fixed distance R is known to follow a power-law M~" 
with a* = 1.39 [l3|. Note that, for comparison of our re- 
sults to Araiijo et al. [3], P{M\R) needs to be divided by 
M, as we grow the IP-cluster starting from the outlet at 
the watershed, which is always the highest of the M sites 
of the cluster. Hence we expect P{M\R) ~ M-("*+i\ 
what is indeed in good agreement with our data (see 
Fig. SMS in T^). We can now relate our exponent a of 
the distribution of areas at fixed distance to a* , which de- 
scribes the size-distribution of IP-clusters, as P{A\R) = 
P{(A) {M)\R) oc {Ay" (M) (X M^"/'^/ cx P{M\R). We 
obtain a = ^{a* -f 1) w 2.266, which is very close to 
what we measure {a = 2.23 ± 0.03). Therefore, we can 
relate our results on uncorrelated landscapes to the sub- 
critical point-to-point invasion percolation process [l3| 
and to the mass distribution of avalanches that occur 
during the IP-cluster growth 

In summary, we were able to show that small and 
localized perturbations can have a large impact on the 
shape of watersheds even at very long distances, hence 
having a non-local effect. The distribution of changes 
P{A) is found to decrease as a power-law with exponent 
/3 = 1.65±0.15 on all studied real landscapes from moun- 
tainous (e.g. Rocky Mountains) to rather flat (e.g. US- 
Canadian border). By applying perturbations to model 
landscapes with long-range correlations, we determined 
the dependence of the scaling exponents on the Hurst 
exponent, finding good quantitative agreement with real 
landscapes, for which 0.3 < H < 0.5. The obtained ex- 
ponents a, P and p are independent of the perturbation 
strength A. For uncorrelated landscapes, we derived a 
relation with invasion percolation. It is known that wa- 
tersheds [1J| on uncorrelated landscapes are related to 
"strands" in Invasion Percolation [2^ , random polymers 
in strongly disordered media [24|, paths on MST's (25j . 
the backbone of the optimal path crack 126| and the clus- 
ter perimeter in explosive percolation [27|. Hence, our 
results can be potentially applied to all these problems. 
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